#!/usr/bin/perl
# 2013/05/08 Only Extract AA or BB combination


use strict;
use warnings;
use IO::File;
use Getopt::Std;
use File::Basename;
use Carp qw(croak);
$| = 0;


if(@ARGV != 4 or $ARGV[0] eq '--help'){
	die "Usage : $0 -i [map file] -o [output filtered map file]\n";
}

my %opt;
getopts("i:o:", \%opt);
my $FILE_map = $opt{i};
my $FILE_out = $opt{o};


my $COUNTER_AAorBB = 0;
my $COUNTER_ABorBA = 0;

# id chr1 loc1 direction1 quality1 hind1 hind_loc1 chr2 loc2 direction2 quality2 hind2 hind_loc2 type(unique or repeat) number(alignable locations) locations
my $fh_out = IO::File->new($FILE_out, 'w') or die "cannot write $FILE_out: $!";
#96241027_ChIA-PET_linkerBB      chr1    69      +       0       NA      NA      chr2    4532782 -       37      NA      NA      R:U     2:1     chr2,-4532454,25M,0;:NA
#23307685_ChIA-PET_linkerAA      chr1    69      +       0       NA      NA      chr2    4533340 -       37      NA      NA      R:U     2:1     chr2,-4532454,25M,0;:NA
#166557348_ChIA-PET_linkerAA     chr1    70      -       0       NA      NA      chr2    2716958 -       37      NA      NA      R:U     2:1     chr2,+4532478,25M,0;:NA
#39628447_ChIA-PET_linkerBB      chr1    70      -       0       NA      NA      chr3    330028  +       37      NA      NA      R:U     2:1     chr2,+4532478,25M,0;:NA
#101751072_ChIA-PET_linkerAA     chr1    74      +       0       NA      NA      chr1    11727   -       0       NA      NA      R:R     2:3     chr2,-4532449,25M,0;:chr2,+4520823,25M,0;chr1,+5575996,25M,0;
#119881221_ChIA-PET_linkerAA     chr1    74      +       0       NA      NA      chr1    5575996 +       0       NA      NA      R:R     2:3     chr2,-4532449,25M,0;:chr1,-11702,25M,0;chr2,+4520823,25M,0;
#87340779_ChIA-PET_linkerAA      chr1    74      +       0       NA      NA      chr2    4520823 +       0       NA      NA      R:R     2:3     chr2,-4532449,25M,0;:chr1,-11702,25M,0;chr1,+5575996,25M,0;
#148627531_ChIA-PET_linkerBB     chr1    74      +       0       NA      NA      chr2    4532614 -       37      NA      NA      R:U     2:1     chr2,-4532449,25M,0;:NA
my $fh_map = IO::File->new($FILE_map) or die "cannot open $FILE_map: $!";
while($_ = $fh_map->getline()){
	s/\r?\n//;
	my ($id, $chr1, $loc1, $direction1, $mapQ1, $hinID1, $hinLoc1, $chr2, $loc2, $direction2, $mapQ2, $hinID2, $hinLoc2, $type, $num, $locations) = split /\t/;

	# both side should have data
	if($chr1 eq 'NA' or $chr2 eq 'NA'){
		next;
	}

	if($id =~ m/linkerAA/ or $id =~ m/linkerBB/){
		$COUNTER_AAorBB++;
		$fh_out->print("$_\n");
	}else{
		$COUNTER_ABorBA++;
	}

}
$fh_map->close();
$fh_out->close();

my $SUM = $COUNTER_AAorBB + $COUNTER_ABorBA;
print "Linker combination\n";
printf "AAorBB: %d (%.3f)\n", $COUNTER_AAorBB, $COUNTER_AAorBB / $SUM * 100;
printf "ABorBA: %d (%.3f)\n\n", $COUNTER_ABorBA, $COUNTER_ABorBA / $SUM * 100;



